{
 "cells": [
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "Text provided under a Creative Commons Attribution license, CC-BY.  All code is made available under the FSF-approved BSD-3 license.  (c) Lorena A. Barba, Gilbert F. Forsyth 2017. Thanks to NSF for support via CAREER award #1149784."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[@LorenaABarba](https://twitter.com/LorenaABarba)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "12 steps to Navier–Stokes\n",
    "======\n",
    "***"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You should have completed Steps [1](./01_Step_1.ipynb) and [2](./02_Step_2.ipynb) before continuing. This Jupyter notebook continues the presentation of the **12 steps to Navier–Stokes**, the practical module taught in the interactive CFD class of [Prof. Lorena Barba](http://lorenabarba.com). "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 3: Diffusion Equation in 1-D\n",
    "-----\n",
    "***"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The one-dimensional diffusion equation is:\n",
    "\n",
    "$$\\frac{\\partial u}{\\partial t}= \\nu \\frac{\\partial^2 u}{\\partial x^2}$$\n",
    "\n",
    "The first thing you should notice is that —unlike the previous two simple equations we have studied— this equation has a second-order derivative. We first need to learn what to do with it!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Discretizing $\\frac{\\partial ^2 u}{\\partial x^2}$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The second-order derivative can be represented geometrically as the line tangent to the curve given by the first derivative.  We will discretize the second-order derivative with a Central Difference scheme: a combination of Forward Difference and Backward Difference of the first derivative.  Consider the Taylor expansion of $u_{i+1}$ and $u_{i-1}$ around $u_i$:\n",
    "\n",
    "$u_{i+1} = u_i + \\Delta x \\frac{\\partial u}{\\partial x}\\bigg|_i + \\frac{\\Delta x^2}{2} \\frac{\\partial ^2 u}{\\partial x^2}\\bigg|_i + \\frac{\\Delta x^3}{3!} \\frac{\\partial ^3 u}{\\partial x^3}\\bigg|_i + O(\\Delta x^4)$\n",
    "\n",
    "$u_{i-1} = u_i - \\Delta x \\frac{\\partial u}{\\partial x}\\bigg|_i + \\frac{\\Delta x^2}{2} \\frac{\\partial ^2 u}{\\partial x^2}\\bigg|_i - \\frac{\\Delta x^3}{3!} \\frac{\\partial ^3 u}{\\partial x^3}\\bigg|_i + O(\\Delta x^4)$\n",
    "\n",
    "If we add these two expansions, you can see that the odd-numbered derivative terms will cancel each other out.  If we neglect any terms of $O(\\Delta x^4)$ or higher (and really, those are very small), then we can rearrange the sum of these two expansions to solve for our second-derivative.  \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$u_{i+1} + u_{i-1} = 2u_i+\\Delta x^2 \\frac{\\partial ^2 u}{\\partial x^2}\\bigg|_i + O(\\Delta x^4)$\n",
    "\n",
    "Then rearrange to solve for $\\frac{\\partial ^2 u}{\\partial x^2}\\bigg|_i$ and the result is:\n",
    "\n",
    "$$\\frac{\\partial ^2 u}{\\partial x^2}=\\frac{u_{i+1}-2u_{i}+u_{i-1}}{\\Delta x^2} + O(\\Delta x^2)$$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Back to Step 3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now write the discretized version of the diffusion equation in 1D:\n",
    "\n",
    "$$\\frac{u_{i}^{n+1}-u_{i}^{n}}{\\Delta t}=\\nu\\frac{u_{i+1}^{n}-2u_{i}^{n}+u_{i-1}^{n}}{\\Delta x^2}$$\n",
    "\n",
    "As before, we notice that once we have an initial condition, the only unknown is $u_{i}^{n+1}$, so we re-arrange the equation solving for our unknown:\n",
    "\n",
    "$$u_{i}^{n+1}=u_{i}^{n}+\\frac{\\nu\\Delta t}{\\Delta x^2}(u_{i+1}^{n}-2u_{i}^{n}+u_{i-1}^{n})$$\n",
    "\n",
    "The above discrete equation allows us to write a program to advance a solution in time. But we need an initial condition. Let's continue using our favorite: the hat function. So, at $t=0$, $u=2$ in the interval $0.5\\le x\\le 1$ and $u=1$ everywhere else. We are ready to number-crunch!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAGdCAYAAADAAnMpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABDr0lEQVR4nO3deXxU5b0/8M+ZPdtMyEpCFsIW9hCEREQFFEW0VNpbBTdc22r1V5Vu8movlGvvRVuXtrfUbir2iiAqYuuCIhBQRAKBsO9JSEJC9mSyTmZ5fn9MZiBCQibM5Jnl83695gUzOZP5npxJzmfO8z3PUYQQAkRERESSqGQXQERERKGNYYSIiIikYhghIiIiqRhGiIiISCqGESIiIpKKYYSIiIikYhghIiIiqRhGiIiISCqN7AL6wuFwoKKiAlFRUVAURXY5RERE1AdCCDQ3NyM5ORkqVc/HPwIijFRUVCA1NVV2GURERNQPZWVlSElJ6fHrARFGoqKiADhXxmg0Sq6GiIiI+sJsNiM1NdW9H+9JQIQR19CM0WhkGCEiIgowl2uxYAMrERERScUwQkRERFIxjBAREZFUDCNEREQkFcMIERERScUwQkRERFIxjBAREZFUDCNEREQkFcMIERERScUwQkRERFIxjBAREZFUDCNEREQkFcMIkR/bV9qAV/JO42B5k+xSiIh8JiCu2ksUSoQQyDtRg7/kncau4noAwPMApo+IxaMzhuPaEXGXvQImEVEgYRgh8hM2uwMfHqjEX7adxrFzzQAArVrBlPQY7C6px45Tddhxqg7jhxjxw+uH49YJSVCrGEqIKPApQgghu4jLMZvNMJlMaGpqgtFolF0OkVe1d9rx9u5S/P2LYpxtbAcAROjUuDs3DQ9dm4EkUxjONrbjH18UYW1+GdqtdgBAWkw4vn/9MNxxVQoMWrXMVSAiuqS+7r8ZRogkaWjtxBs7S/DGVyVoaLMCAOIidXhwegbuzU2HKVzb5+c8cM1Q3Hf10Es+h4hIFoYRIj9V3dyBP289jbd3dz/K8YPrh+F7fTzK0dZpw7rdZRcdTbkrJw2PzRyO2Ei9T9eBiKgvGEaI/FCLxYbb/vgFztS1AQDGJRvx2MzhmDu+f/0fVrsDH32jz2REQiQ+/H/XcuiGiKTr6/6bDaxEA+jZfx/Bmbo2JJsMeP57E6/4zBitWoX52UNw+6Rk5J2owc/fPYBT1S147pNj+PW3x3mxciIi3+E8I0QD5NPD5/D2njIoCvDSgkm4bmS8107RVRQFszIT8MIdWQCAVV+VYPuJGq98byIiX2MYIRoA1c0dWLL+IADgB9cNw9XDYn3yOjNGxeP+aekAgJ++sx8NrZ0+eR0iIm9iGCHyMSEEfvHuAdS3dmJMkhGLbx7l09d7Zu4YDI+PQHWzBb/ccBAB0BZGRCGOYYTIx1bvKsXW4zXQaVT4/YJJ0Gt821gaplPj9wuyoVEp+PjgOazfe9anr0dEdKUYRoh8qKimBf/90VEAwM/nZCJzcNSAvO6EFBOemj0SALDsX4dRVt82IK9LRNQfDCNEPmK1O/D024Vot9oxfUQsHpqeMaCv/+iM4bgqfRBaLDb8ZN1+2B0criEi/8QwQuQjf9pyCvvLm2A0aPDCHVlQDfB1ZDRqFV6+cxIidGrkl9Tjb9uLBvT1iYj6imGEyAf2lTbgT1tPAQB+850JSDKFSakjLTYcy+Y55xt5adNxHDrbJKUOIqLeMIwQeVmrxYan3y6E3SFw+6RkfDsrWWo9d0xJwc1jE2G1Czz9diE6uqagJyLyFwwjRF72m4+OoqSuDUkmA/7r2+NllwNFUbDiuxMQF6nHyeoWPL/xmOySiIi6YRgh8qLPj1RhTX4pAODFO7L85iq6sZF6/O57EwEAr+8owRcnOTsrEfkPhhEiL6ltseCZ9QcAAI9cm4FrRsRJrqi7WaMTcO/VaQCcs7M2tnF2ViLyDwwjRF4ghMAz7x1AbUsnRg+Owk/nZMou6ZJ+eetYDIuLQJXZgl9uOMTZWYnILzCMEHnB2t1l+PxoNXRqFV5eMAkGrW9nWe2vMJ0aLy+YBI1KwUcHKrGhkLOzEpF8DCNEV6i+tRPPfngEAPDTOaMwJskouaLeZaVG48c3OmdnXfrBYZg7rJIrIqJQxzBCdIXe2VOGtk47xiYZ8ci1w2SX0yc/mjkcIxIi0dxhw4Z9PDpCRHIxjBBdAYdD4K2us2cWTUsf8FlW+0ujVuHeXGcz65tfn2HvCBFJxTBCdAW+PFWLM3VtiNJr8O1Jcic389R3r0pBmFaNE1Ut2HOmQXY5RBTCGEaIrsCbX58BAHx38hCE6zSSq/GM0aB1zw7rWg8iIhkYRoj6qbKpHZuPVQMA7rk6XXI1/XNvV92fHDyHuhaL5GqIKFQxjBD109r8MtgdAjkZMRiVGCW7nH6ZkGLCxBQTOu0OvFNQLrscIgpRDCNE/WCzO7B2t7Nx9Z6uRtBAdW+u8+jIW7tK4XCwkZWIBh7DCFE/fH60GlVmC2IjdLhl/GDZ5VyRb2UlIcqgQWl9G744VSu7HCIKQQwjRP2wepez4fOOKanQa/xzttW+Ctdp8B+TUwAAq9nISkQSMIwQeaikthVfnKyFogB35wT2EI2La6jp86NVqGxql1wNEYUahhEiD63pmuTs+pHxSIsNl1yNd4xMjEJuRgwcwtmYS0Q0kBhGiDzQYbVj3R7nzvreAD2dtyeu05PX7i6F1e6QXA0RhRKGESIPfHKoEg1tViSZDJiVGS+7HK+6ZdxgxEXqUGW2YPPRKtnlEFEIYRgh8sDqr51DNHflpEGjDq5fH51GhTunpAIAVu8qlVwNEYWS4PprSuRDx86ZsedMA9QqBQunpsouxyfuykmDogBfnKxFcW2r7HKIKEQwjBD1keuoyM1jE5FgNEiuxjdSY8Ixc5Rz+MnVqEtE5GsMI0R90Gqx4f19ZwEEX+PqN93TNSPrO3vK0GG1S66GiEIBwwhRH3xQWIEWiw0ZcRGYNixWdjk+NWt0ApJNBjS0WfHJoUrZ5RBRCGAYIboMIQTe7JqZ9J7cNKhUiuSKfEutUnBX12Rub37NoRoi8j2GEaLLKCxrxJFKM3QalXva9GC3YGoqNCoFBWcacLTSLLscIgpyDCNEl+E6OvCtiUkYFKGTXM3ASDAacPO4RADnr8NDROQrDCNEvWhs68SHByoAnG/sDBX3dq3v+3vPosVik1wNEQUzhhGiXrxbUA6LzYExSUZMTouWXc6AmjY8FsPiItDaaccHhWdll0NEQYxhhKgHQgi81TUT6b1Xp0FRgrtx9ZsURcHduecbWYUQkisiomDFMELUg52n61BU24oInRq3TxoiuxwpvndVCvQaFY5WmrGvrFF2OUQUpBhGiHrguj7LdyYPQaReI7kaOaLDdfjWxGQA52egJSLyNoYRokuoNnfg08PnAIRe4+o33XO1c6jmwwMVaGzrlFwNEQUjj8PI9u3bMW/ePCQnJ0NRFGzYsOGyz1m9ejWysrIQHh6OpKQkPPTQQ6irq+tPvUQD4t8HKmFzCExOi8aYJKPscqTKTnX+DCw2Bz4+eE52OUQUhDwOI62trcjKysLKlSv7tPyOHTuwaNEiPPzwwzh8+DDeeecd5Ofn4/vf/77HxRINlK3HqgEAt05IklyJfIqi4NbxgwEAW7p+LkRE3uTxQPjcuXMxd+7cPi+/c+dODB06FD/+8Y8BABkZGfjhD3+I559/3tOXJhoQLRYbdhU7j9zNGp0guRr/MGt0Al7cdAI7TtWiw2qHQauWXRIRBRGf94xMmzYNZWVl+PjjjyGEQFVVFd59913ceuutPT7HYrHAbDZ3uxENlC9P1sJqF0iPDcewuAjZ5fiFcclGJBr1aLfasau4XnY5RBRkfB5Gpk+fjtWrV2PBggXQ6XQYPHgwTCZTr8M8K1asgMlkct9SU1N9XSaRW95x51DErMyEkJtbpCeKomBWpvMo0VYO1RCRl/k8jBw5cgRPPvkkli5dioKCAmzcuBElJSV49NFHe3zOkiVL0NTU5L6VlZX5ukwiAM6JzrZ2hZEbOETTjWvIauvxak6ARkRe5fPJE1asWIHp06fjZz/7GQBg4sSJiIiIwHXXXYff/OY3SEq6uEFQr9dDr9f7ujSiixyuMKPKbEGYVo3cYTGyy/Er146Ig1at4ExdG4pqWzE8PlJ2SUQUJHx+ZKStrQ0qVfeXUaudzW/8dEX+xjUEMX1EHPQaNmleKEKvQW5GLAAO1RCRd3kcRlpaWlBYWIjCwkIAQHFxMQoLC1Fa6pydccmSJVi0aJF7+Xnz5mH9+vV45ZVXUFRUhB07duDHP/4xcnJykJyc7J21IPKSLRyi6ZVrqIan+BKRN3kcRvbs2YPs7GxkZ2cDABYvXozs7GwsXboUAFBZWekOJgDwwAMP4KWXXsKf/vQnjB8/HnfccQcyMzOxfv16L60CkXfUt3aisOv6K7NGx8stxk+5Qlp+cT2aO6ySqyGiYKGIABgrMZvNMJlMaGpqgtEY2rNhku+8v68cT7+9H2OSjPjkyetkl+O3Zr2Qh+LaVrxyz2TM5aRwRNSLvu6/eW0aoi5bjtUAAG7gUZFeuU/xPc6hGiLyDoYRIgA2uwPb2C/SJze4T/GtgcPh9wdWiSgAMIwQAdhb2ghzhw3R4VpMSh0kuxy/NjVjECJ0atQ0W3C4grMjE9GVYxghwvmzQ2aMiodaxVlXe6PXqDF9RBwAnlVDRN7BMEKE81PAc4imb1w/py3sGyEiL2AYoZB3trEdx841Q6U4j4zQ5bnmGzlQ3ojaFovkaogo0DGMUMhzzSY6OW0QosN1kqsJDIlGA8YlGyEEkHe8RnY5RBTgGEYo5LnCyCwO0XjkhtE8xZeIvINhhEJah9WOHadrAbBfxFOu8Lb9RA2sdofkaogokDGMUEjbWVSHDqsDSSYDRg+Okl1OQMlKiUZMhA7NHTYUnGmQXQ4RBTCGEQppeV1DNDMzE6AoPKXXE2qV4m745VV8iehKMIxQyBJC8Cq9V4hX8SUib2AYoZB1uqYFZfXt0GlUmD4iVnY5AWnGSOckcSerW1BW3ya7HCIKUAwjFLJcn+avHhaLcJ1GcjWByRSuxVVpzunzeVYNEfUXwwiFLFcYuSGTE51dCddQDftGiKi/GEYoJJk7rNhT4jwDhPOLXJlZo51h7qvTdWjvtEuuhogCEcMIhaQvT9bC5hAYFh+B9NgI2eUEtMzEKCSbDLDYHNhZVCu7HCIKQAwjFJLOD9HwqMiVUhSFZ9UQ0RVhGKGQ43AIXqXXy9xTwx+rgRBCcjVEFGgYRijkHDzbhNqWTkTqNZgyNEZ2OUHhmuFx0GtUONvYjhNVLbLLIaIAwzBCIcc1lHDdyDjoNPwV8IYwnRrThjvnauFQDRF5in+JKeS45sPgWTTexav4ElF/MYxQSKlptuBAeRMAYCbnF/GqWV3NwAVnGtDUZpVcDREFEoYRCimuxtUJQ0xIiDJIria4pMaEY0RCJOwOge0na2SXQ0QBhGGEQgqHaHzrBs7GSkT9wDBCIcNqd+CLE85JuXhKr2+4hmryTtTA7uApvkTUNwwjFDJ2l9Sj2WJDXKQOE4eYZJcTlKYMHYQogwb1rZ3YX94ouxwiChAMIxQy8o47+xhmjEqASqVIriY4adUqXD/S2RjMoRoi6iuGEQoZO0/XAQCuHxUnuZLg5vr5fl1UJ7kSIgoUDCMUEpo7rDhc4TylNzcjVnI1wc31891f1oQOK6/iS0SXxzBCIaHgTAMcAkiPDcdgE0/p9aX02HAkROnRaXegsKxRdjlEFAAYRigk5BfXAwCm8lo0PqcoCqZmOH/Orp87EVFvGEYoJOwuce4UczIYRgZCbtfP2fVzJyLqDcMIBb0Oqx37y1z9IgwjA8EV+grONMBqd0iuhoj8HcMIBb3CskZ02h1INOqRFhMuu5yQMCohCqYwLdo67ThcYZZdDhH5OYYRCnquvoWcjFgoCucXGQgqleLuz8kv5im+RNQ7hhEKeu4wMnSQ5EpCS06G8+fNJlYiuhyGEQpqVrsDBWcaADiPjNDAcf28d5c0wMHr1BBRLxhGKKgdrjCj3WpHdLgWIxMiZZcTUsYlGxGuU6Op3YoT1c2yyyEiP8YwQkHN1a8wdWgMr0czwLRqFa5K51ANEV0ewwgFNddOkKf0ypHT1cS6i2GEiHrBMEJBy+EQnHlVsgtnYhWCfSNEdGkMIxS0jlc1w9xhQ7hOjXHJRtnlhKRJqdHQqVWoabagpK5NdjlE5KcYRihouaYivyp9EDRqvtVlMGjVyEo1AQB2c6iGiHrAv9AUtHaxX8QvuKaGZ98IEfWEYYSCkhCi28yrJI/r559fwplYiejSGEYoKJXUtaGm2QKdWoWJKSbZ5YS0yWnRUClAWX07KhrbZZdDRH6IYYSCkmt+kUmp0TBo1ZKrCW1RBi3GJXf1jZRwqIaILsYwQkFpl3uIhv0i/oB9I0TUG4YRCkquT+AMI/7BtR14Rg0RXQrDCAWdisZ2lNW3Q61SMDmdV+r1B65J505Wt6CuxSK5GiLyNwwjFHRcR0XGJxsRqddIroYAICZCh1GJzgsV7i5pkFwNEfkbhhEKOrs4Bbxfcm0PXjSPiL6JYYSCTj6bV/2Sa3twvhEi+iaGEQoqdS0WnKpuAcAjI/7GFUaOVJjR3GGVXA0R+ROGEQoqrn6EzMQoDIrQSa6GLpRkCkNaTDgcAig4w74RIjqPYYSCCodo/Jt7qIZ9I0R0AYYRCiqufoSpDCN+KYdNrER0CQwjFDTMHVYcqTADOL/TI//iOjKyv7wRHVa75GqIyF8wjFDQKDjTAIcA0mPDMdhkkF0OXUJ6bDgSovSw2gX2lTbKLoeI/ITHYWT79u2YN28ekpOToSgKNmzYcNnnWCwW/PKXv0R6ejr0ej2GDh2K1157rT/1EvXINdU4j4r4L0VRzk8Nz4vmEVEXj6enbG1tRVZWFh566CF897vf7dNz7rzzTlRVVeHVV1/FiBEjUFlZCYfD4XGxRL1h82pgyM2IwYcHKtk3QkRuHoeRuXPnYu7cuX1efuPGjdi2bRuKiooQE+PcSQwdOtTTlyXqVYfVjv3ljQAYRvydq7m44EwDrHYHtGqOFhOFOp//FfjXv/6FKVOm4Le//S2GDBmCUaNG4ac//Sna29t7fI7FYoHZbO52I+rNvtJGWO0CiUY90mLCZZdDvRiVEAVTmBbtVjsOnW2SXQ4R+QGfh5GioiJ8+eWXOHToEN5//338/ve/x7vvvosf/ehHPT5nxYoVMJlM7ltqaqqvy6QAd36IJhaKokiuhnqjUim8Tg0RdePzMOJwOKAoClavXo2cnBzceuuteOmll/DGG2/0eHRkyZIlaGpqct/Kysp8XSYFOFczJIdoAkMum1iJ6AI+v756UlIShgwZApPJ5H5szJgxEEKgvLwcI0eOvOg5er0eer3e16VRkLDaHe7pxXMZRgLChTOxOhwCKhWPZhGFMp8fGZk+fToqKirQ0tLifuzEiRNQqVRISUnx9ctTCDh0tgntVjuiw7UYER8puxzqg3HJRoTr1DB32HC8qll2OUQkmcdhpKWlBYWFhSgsLAQAFBcXo7CwEKWlpQCcQyyLFi1yL3/33XcjNjYWDz74II4cOYLt27fjZz/7GR566CGEhYV5Zy0opLn6DqYOjeEn7AChUatwVfogAOwbIaJ+hJE9e/YgOzsb2dnZAIDFixcjOzsbS5cuBQBUVla6gwkAREZGYtOmTWhsbMSUKVNwzz33YN68efjjH//opVWgUOfamXGIJrDwOjVE5OJxz8jMmTMhhOjx66tWrbrosdGjR2PTpk2evhTRZTkcgs2rAcrdN1JSDyEEz4IiCmGcbYgC2vGqZpg7bIjQqTE2ySi7HPJAVmo0dGoVapotKKlrk10OEUnEMEIBzXWI/6qhMdBwJs+AYtCqMSk1GgCQX1wntxgikop/vSmguSc7GzpIciXUH1MznNttF/tGiEIawwgFLCEE8kvOn0lDgScnIxYAm1iJQh3DCAWssvp21DRboFUryOo63E+BZXJaNBQFKG9oR5W5Q3Y5RCQJwwgFrIJS56fp8UNMMGjVkquh/ogyaJGZGAUA2Ns1iy4RhR6GEQpYringr0pjv0ggc01+VsAwQhSyGEYoYBWcaQRwfmdGgckdRkoZRohCFcMIBaTmDiuOnzMDACYzjAQ0Vxg5dLYJHVa75GqISAaGEQpI+8ua4BBAyqAwJBoNssuhK5AWE464SB2sdoFDZ5tkl0NEEjCMUEBy94vwqEjAUxQFk9PYN0IUyhhGKCC5+gsYRoIDm1iJQhvDCAUch0NgX9dOazLPpAkKrjCyt7Sh1wtxElFwYhihgHOyugXNFhvCdWqMHhwluxzygvFDTNCqFdS2dKK0nhfNIwo1DCMUcFyH8ielRvPieEHCoFVj/BATAA7VEIUi/iWngMPm1eB0FZtYiUIWwwgFnL1dzaucXyS4sImVKHQxjFBAqWuxoLi2FQAwOZVhJJi4wuXxqmY0d1glV0NEA4lhhALK3tJGAMDIhEiYwrVyiyGvSjQakDIoDEIAhWWNssshogHEMEIBhf0iwY1DNUShiWGEAorrMvPsFwlODCNEoYlhhAJGp82B/eWNAHhkJFi5JrErLG2E3cHJz4hCBcMIBYwjlWZYbA5Eh2sxLC5CdjnkA6MHRyFcp0azxYaT1c2yyyGiAcIwQgHD3S+SNgiKokiuhnxBo1ZhUmo0AA7VEIUShhEKGOwXCQ3sGyEKPQwjFBCEENhzph4A+0WCnSts7mUYIQoZDCMUECqaOlBltkCtUpCVEi27HPIh12R2JXVtqG2xSK6GiAYCwwgFBNch+3HJRoTp1JKrIV8yhWsxMiESAI+OEIUKhhEKCO5+kTQO0YQCd99IKcMIUShgGKGAwJlXQwv7RohCC8MI+b22ThuOVJoBMIyEiild23l/eRM6bQ7J1RCRrzGMkN/bX9YEu0MgyWRAcnSY7HJoAGTERWBQuBadNgcOVzTJLoeIfIxhhPze3lLOLxJqFEXhfCNEIYRhhPzehTOvUuhw942wiZUo6DGMkF9zOIR7Z8R+kdDiCp8FZxogBC+aRxTMGEbIrxXVtqKxzQqDVoWxyUbZ5dAAmpgSDY1KQZXZgrON7bLLISIfYhghv+Y6tXNiSjS0ar5dQ0mYTo1xXQGUfSNEwY1/3cmvcX6R0Mb5RohCA8MI+TXXDJxsXg1NnImVKDQwjJDfamzrxKnqFgA8rTdUucLI0cpmtFpskqshIl9hGCG/ta+0EQAwLC4CMRE6ucWQFEmmMCSbDLA7BPaXN8ouh4h8hGGE/JarX4RHRUIb+0aIgh/DCPktNq8SAM7EShQCGEbIL1ntDhSWNQJgGAl1V7lnYm2Ew8HJz4iCEcMI+aVjlc1ot9oRZdBgRHyk7HJIojFJRhi0KjS1W1FU2yK7HCLyAYYR8ksFZ+oBAJPTBkGlUiRXQzJp1SpkpUQD4FANUbBiGCG/VNB1Jg2HaAhg3whRsGMYIb+0l82rdAGGEaLgxjBCfqeyqR1nG9uhUoCs1GjZ5ZAfyO6agfd0TSsaWjslV0NE3sYwQn5n75lGAMDowUZE6jVyiyG/EBOhw7D4CADAvjIeHSEKNgwj5Hc4vwhdiuv6RByqIQo+DCPkd9wXx2MYoQuwb4QoeDGMkF/psNpx+GwTAIYR6s71fthf1gSr3SG5GiLyJoYR8iv7yxphcwjER+mRMihMdjnkR4bHR8Jo0KDdaseRCrPscojIixhGyK/kFzsnO8sZGgNF4WRndJ5KpWDq0BgAwO6SesnVEJE3MYyQX8nv2snkZMRIroT8ket94QqtRBQcGEbIb9jsDndzousTMNGFpmacPzLCi+YRBQ+GEfIbhyvMaOu0w2jQIHNwlOxyyA+NTzYhTKtGQ5sVp2p40TyiYMEwQn7D1QcwdWgM1Lw4Hl2CTqPC5PRoAByqIQomDCPkN3YVs1+ELi9naCwAhhGiYOJxGNm+fTvmzZuH5ORkKIqCDRs29Pm5O3bsgEajwaRJkzx9WQpyDoc4f2SEYYR6MTXDOd9IfnE9hGDfCFEw8DiMtLa2IisrCytXrvToeY2NjVi0aBFuvPFGT1+SQsDJ6hY0tlkRplVjfLJJdjnkx7JTB0GrVnDO3IGy+nbZ5RCRF3h8FbK5c+di7ty5Hr/Qo48+irvvvhtqtdqjoykUGlyn9E5Oj4ZOw9FD6lmYTo2JKdEoONOA/JJ6pMWGyy6JiK7QgPzVf/3111FUVIRly5b1aXmLxQKz2dztRsHNNf7PU3qpL1zvk/ziOsmVEJE3+DyMnDx5Es888wzefPNNaDR9OxCzYsUKmEwm9y01NdXHVZJMQgj3ToXNq9QXuZz8jCio+DSM2O123H333Vi+fDlGjRrV5+ctWbIETU1N7ltZWZkPqyTZyurbUWW2QKtWkJ3Ki+PR5V01dBAUBSipa0O1uUN2OUR0hTzuGfFEc3Mz9uzZg3379uGJJ54AADgcDgghoNFo8Nlnn+GGG2646Hl6vR56vd6XpZEf2dV1VGTCEBPCdGrJ1VAgMBq0GDPYiCOVZuSX1ONbE5Nll0REV8CnYcRoNOLgwYPdHvvzn/+MLVu24N1330VGRoYvX54ChPvieBmxkiuhQJKTEeMMI8UMI0SBzuMw0tLSglOnTrnvFxcXo7CwEDExMUhLS8OSJUtw9uxZ/POf/4RKpcL48eO7PT8hIQEGg+Gixyl0ueYXyWW/CHkgNyMGq74qYd8IURDwOIzs2bMHs2bNct9fvHgxAOD+++/HqlWrUFlZidLSUu9VSEGtytyBkro2KIqzD4Cor1yT4x2vakZjWyeiw3WSKyKi/lJEAExhaDabYTKZ0NTUBKPRKLsc8qJ/76/A/1uzD2OTjPj4yetkl0MB5oYX81BU04p/LJqC2WMTZZdDRN/Q1/03Z5ciqfJ5PRq6Au5TfEs4VEMUyBhGSCr2i9CVyOF8I0RBgWGEpGls68Sxc80AgCmceZX6wTUT66GzTWi12CRXQ0T9xTBC0uwuaQAADIuPQHwU55Uhz6UMCseQ6DDYHAL7Shtll0NE/cQwQtJwiIa8IYd9I0QBj2GEpNnF5lXygvN9I7xoHlGgYhghKVotNhw62wSAV+qlK+N6/+wrbYTFZpdcDRH1B8MISbGvtBF2h8CQ6DCkDAqXXQ4FsOHxEYiN0MFic7gDLhEFFoYRksJ1SJ1DNHSlFEVxv4928RRfooDEMEJSuHYaHKIhb3C9jzjfCFFgYhihAWex2bGvrBEAj4yQd7jeRwUlDbA7/P4KF0T0DQwjNOAOljeh0+ZAbIQOw+MjZJdDQWBMkhFReg2aLTYcrTTLLoeIPMQwQgPuwlN6FUWRXA0FA7VKwZSuqz5zqIYo8DCM0IDLZ78I+cBUXqeGKGAxjNCAsjsECs44p4Fnvwh5k2sm390l9RCCfSNEgYRhhAbU0UozWiw2ROk1GJNklF0OBZEJQ6Kh16hQ19qJ0zWtssshIg8wjNCAcvWLXDV0ENQq9ouQ9+g0KmSnRQPgUA1RoGEYoQG1m9ejIR/KyYgFcP4ijEQUGBhGaMAIIdxXVuWVeskXctnEShSQGEZowJyuaUF9ayf0GhUmDImWXQ4Foey0aGhUCs42tqO8oU12OUTURwwjNGBc/SLZadHQafjWI+8L12kwfogJAI+OEAUS7hFowJzvF4mVXAkFswtP8SWiwMAwQgNCCOE+MsJ+EfIlXsGXKPAwjNCAKG9oR2VTBzQqxX36JZEvTEmPgaIARTWtqGm2yC6HiPqAYYQGhOuQ+fghJoTrNJKroWBmCtciMzEKALCHQzVEAYFhhAZEPodoaADlcqiGKKAwjNCA4MXxaCDxonlEgYVhhHyuptmCotpWKArDCA2MnK732dFzZpg7rJKrIaLLYRghn3P1i2QmRsEUrpVcDYWCBKMBGXEREIJ9I0SBgGGEfO6Lk7UAgKuHcX4RGjiuvpEvT9ZJroSILodhhHxKCIG849UAgBmZ8ZKroVAys+v95nr/EZH/Yhghnzpe1YzKpg4YtCpM45ERGkDTR8RBq1ZQVNuKktpW2eUQUS8YRsinthxzfiq9ZngcDFq15GoolEQZtO6G6a08OkLk1xhGyKfyjtUAAGZxiIYkmJWZAADYerxGciVE1BuGEfKZpjYrCkobAAAzu3YKRANp1mhnCP66qA5tnTbJ1RBRTxhGyGe+OFUDu0NgZEIkUmPCZZdDIWh4fCRSBoWh0+bAztM8q4bIXzGMkM+4+kVmjeZREZJDURTc0PX+c70ficj/MIyQTzgcAtu6xulnsl+EJHL1jeQdr4EQQnI1RHQpDCPkEwfPNqGutROReg2mpHMKeJLn6mGx0GtUONvYjpPVLbLLIaJLYBghn3CdSnndyDjoNHybkTxhOjWmDXfOccOhGiL/xL0E+cRWV78Iz6IhP+DqG9nKMELklxhGyOtqmi3YX94EgFPAk3+YOcoZRvacaeBVfIn8EMMIed32E87G1XHJRiQaDZKrIQLSYsMxPD4CdofAl10XbiQi/8EwQl7n6he5gaf0kh9xDRmyb4TI/zCMkFfZ7A73kRHOukr+xBWO847XwOHgKb5E/oRhhLxqb2kjzB02DArXYlJqtOxyiNymDI1BhE6N2hYLDleYZZdDRBdgGCGvcg3RzBgVD7VKkVwN0Xk6jQrXjowDwKv4EvkbhhHyqq2cAp78GPtGiPwTwwh5TUVjO46da4aiANeP5Cm95H9cIXl/eSPqWiySqyEiF4YR8pq8rmvRZKdGY1CETnI1RBdLNBowNskIIYDtJ2tkl0NEXRhGyGt4Si8FglmjnUftthxjGCHyFwwj5BUWmx07Tjknk+IpveTPXH0j20/UwGZ3SK6GiACGEfKS/OJ6tHXakRClx7hko+xyiHqUnTYI0eFaNLVbUVjWKLscIgLDCHnJ1q5D3rMyE6AoPKWX/JdapbgbrHmKL5F/YBghr8g77jqll2fRkP9j3wiRf2EYoStWUtuKotpWaFQKpo+Ik10O0WVdPzIeigIcrTTjXFOH7HKIQh7DCF0x11GRqUNjEGXQSq6G6PJiI/XuyxXkcaiGSDqGEbpiW7rmF+EpvRRIXGfVsG+ESD6GEboibZ02fF1UB4D9IhRYXGHky5O1sNjskqshCm0MI3RFdp6uQ6fNgZRBYRgeHym7HKI+G5dsRHyUHq2dduwpaZBdDlFI8ziMbN++HfPmzUNycjIURcGGDRt6XX79+vW46aabEB8fD6PRiGnTpuHTTz/tb73kZy6cdZWn9FIgUakUzBzVdYovL5xHJJXHYaS1tRVZWVlYuXJln5bfvn07brrpJnz88ccoKCjArFmzMG/ePOzbt8/jYsm/CCG6zS9CFGhcF87bwr4RIqk0nj5h7ty5mDt3bp+X//3vf9/t/v/8z//ggw8+wL///W9kZ2d7+vLkR05Wt+BsYzv0GhWuHhYruxwij107Mg4alYKimlacqWtFemyE7JKIQtKA94w4HA40NzcjJiamx2UsFgvMZnO3G/kf16HtacNjEaZTS66GyHNGgxZThg4CcP6q00Q08AY8jLzwwgtoaWnBnXfe2eMyK1asgMlkct9SU1MHsELqK16ll4IBT/Elkm9Aw8hbb72F5cuXY926dUhI6HkHtmTJEjQ1NblvZWVlA1gl9YW5w+o+A2HmKIYRClyuvpGdp+vQ3slTfIlkGLAwsnbtWjzyyCNYt24dZs+e3euyer0eRqOx2438y5cna2FzCAyPj0BabLjscoj6bWRCJIZEh8Fic2BnUa3scohC0oCEkTVr1uDBBx/EmjVrcNtttw3ES5KPbTnGIRoKDoqiuCfs23yUQzVEMngcRlpaWlBYWIjCwkIAQHFxMQoLC1FaWgrAOcSyaNEi9/JvvfUWFi1ahBdffBG5ubk4d+4czp07h6amJu+sAQ24VosNGw+dAwDMHpMouRqiK+d6H390sBIdVg7VEA00j8PInj17kJ2d7T4td/HixcjOzsbSpUsBAJWVle5gAgB/+9vfYLPZ8PjjjyMpKcl9e/LJJ720CjTQPjxQgRaLDRlxEcjJ6PmsKKJAcd3IeCSbDGhss+LTw+dkl0MUcjyeZ2TmzJkQQvT49VWrVnW7n5eX5+lLkJ9bk+9sKF4wNZWzrlJQUKsU3DElFX/YfBJr8ktx+6QhsksiCim8Ng155Ng5MwrLGqFVK/jeVSmyyyHymjunpkKlAF8X1aO4tlV2OUQhhWGEPLK266jITWMTERepl1wNkfcMiQ7DjK5r1azdXXqZpYnImxhGqM86rHas31sOAFg4NU1yNUTetzDH+b5+d085Om0OydUQhQ6GEeqzjw9WwtxhQ8qgMFw7Ik52OURed8PoBMRH6VHX2onPj1bJLocoZDCMUJ+5hmgWTk2FSsXGVQo+WrUKd05x9kKtyedQDdFAYRihPjlV3YL8knr3WQdEwWrBFOdQzZenalFW3ya5GqLQwDBCfbK261PirMwEJBoNkqsh8p202HBcOyIOQgBv7+Z1sYgGAsMIXZbFZsd7XY2rd+XwqAgFv4Vd7/N3Cspgs7ORlcjXGEbosj47XIWGNisGGw3uUx+JgtlNYxMRE6FDldmCrcdrZJdDFPQYRuiyXHMu3Dk1FRo13zIU/PQatXtSv7VsZCXyOe5ZqFdn6lqx41QdFAXuswyIQsGCqc6hmq3Hq1HZ1C65GqLgxjBCvVrb1cB3/ch4pAwKl1wN0cAZHh+JnIwYOASwbne57HKIghrDCPXIanfgnT1sXKXQ5Xrfr9tTBruj5wuEEtGVYRihHm0+Wo3aFgviIvW4cUyi7HKIBtzc8UkwhWlxtrEdX5xkIyuRrzCMUI9cM1DeMSUFWjauUggyaNX4TvYQAJyRlciXuIehSypvaMP2rk+CC6dyiIZC111dF8/bfLQa1c0dkqshCk4MI3RJ6/aUQwjgmuGxSI+NkF0OkTSZg6OQnRYNm0Pg3QI2shL5AsMIXcRmd+CdPc6zaFyfColCmev34O3dZXCwkZXI6xhG6CLbTtSgsqkDg8K1uHkcG1eJvjUxCVF6Dc7UtWFnUZ3scoiCDsMIXWRNvvOoyH9MToFeo5ZcDZF84ToNvj0pGQAbWYl8gWGEuqkyd2Dr8WoA5y8WRkTnh2o+O1yF+tZOydUQBReGEermna7JnXKGxmBEQpTscoj8xvghJkwYYkKn3YH1e9nISuRNDCPk5nAI9/TvPCpCdDHX78Vb+aUQgo2sRN7CMEJuX56qRXlDO4wGDW6dkCS7HCK/8+2sZIRp1SiqacXukgbZ5RAFDYYRclu729mY953sITBo2bhK9E1RBi3mZTmD+lo2shJ5DcMIAQBOVDXjs8NVAICFnFuEqEeuRtYPD1SiuLZVcjVEwYFhhCCEwK/ePwSbQ+CmsYkYk2SUXRKR35qUGo3rR8Wj0+7A0g8OsXeEyAsYRgjv7T2L/JJ6hGnVWDZvrOxyiPyaoij4r2+Pg06jwhcna/HRwUrZJREFPIaRENfY1okVHx8FAPz4xpFIGRQuuSIi/zc0LgKPzRgOAPivfx9Bc4dVckVEgY1hJMT97tPjqGvtxMiESDx8bYbscogCxmMzhyM9NhzVzRb8/vOTssshCmgMIyGssKwRb3WdEfDs/PHQafh2IOorg1aN5d8eBwBY9VUJjlSYJVdEFLi49wlRdofAL98/CCGA704egquHxcouiSjgzMxMwK0TBsPuEPjVhoO8oi9RPzGMhKj/21mCwxVmGA0aLJk7RnY5RAFr6bfGIUKnxt7SRqzbUya7HKKAxDASgqrNHXjxsxMAgJ/dMhrxUXrJFREFrsEmA56+aRQA4LmNx3gRPaJ+YBgJQf/98VE0W2zISjHhbk5wRnTF7r9mKEYPjkJjmxXPf3JMdjlEAYdhJMTsOFWLDworoFKA38yfALVKkV0SUcDTqlX4zfzxAIC395RhT0m95IqIAgvDSAix2Oz4zw8OAQDuuzodE1JMkisiCh5ThsbgzikpAIBfbTgEm90huSKiwMEwEkL+vr0IRTWtiIvUY/HNmbLLIQo6z8wdg+hwLY6da8aqr0pkl0MUMBhGQkRZfRv+d8spAMCvbhsDU5hWckVEwScmQodnbhkNAHh50wlUNrVLrogoMDCMhAAhBJb96zAsNgemDYvF7ZOSZZdEFLTunJKK7LRotHba8ZsPj8ouhyggMIyEgM+OVGHLsWpo1QqenT8eisKmVSJfUakU/Gb+eKgU4KODldh2okZ2SUR+j2EkyLV12rD8X4cBAD+4fhhGJERKrogo+I1LNuGBa5zXelr6wSF0WO2SKyLybwwjQe4Pm0+ioqkDQ6LD8MSskbLLIQoZT980EglRepypa8Mreadll0Pk1xhGgthXp2rx6hfFAIDl3x6HMJ1ackVEoSPKoMV/fmssAOCVvNPIL+bcI0Q9YRgJUp8ePocHXt8Nm0Pg1gmDMXtsouySiELOtyYm4cbRCei0O7DotV3YerxadklEfolhJAi9V1COH63ei067A3PGJeLlBZNkl0QUkhRFwZ/unoxZmfHosDrw/Tf24N/7K2SXReR3GEaCzKodxfjJO/thdwh876oUrLx7MvQaDs8QyRKmU+Ov903BvKxk2BwCP167D2vyS2WXReRXGEaChBACf/j8JH797yMAgIemZ+C3/zERGjU3MZFsOo0Kv18wCffkpkEIYMn6g/jLNja1ErloZBdAV87hEPjNR0fx2g5ns+rTs0fhxzeO4HwiRH5E3TX/iClMiz/nncZznxxDU7sVP5+Tyd9VCnkMIwHOZnfgmfUH8W5BOQBg2byxeHB6huSqiOhSFEXBz28ZDWOYFs99cgyv5J2Gud2KZ28fDxWvoE0hjMfwA5jFZscTb+3DuwXlUKsUvHhHFoMIUQB4dMZw/M93JkBRgNW7SvHU24Ww8iq/FMJ4ZCRAtVpsePTNAnxxshY6tQr/e3c25owbLLssIuqju3PTEGXQ4Om3C/Gv/RVosdjw53smw6BlwzmFHh4ZCUBNbVbc++oufHGyFuE6NV5/cCqDCFEAmpeVjL/fPwUGrQpbjlVj0Wv5aO6wyi6LaMAxjASYopoWLPjbTuwrbYQpTIvVj+Ri+og42WURUT/NykzAPx/KRZReg/zietz1969RWtcmuyyiAaUIIYTsIi7HbDbDZDKhqakJRqNRdjlSHCxvwl+2ncbHhyohBJAQpcf/PZyLzMFRsksjIi84dLYJ97+Wj7rWTqgU4FsTk/HojOEYmxyaf/MoOPR1/80w4seEEPjqdB1eyTuNL0/Vuh+/cXQCfv3tcUiNCZdYHRF5W3FtK379r8PYdqLG/djMzHg8NmM4cjJieAowBRyGkQBmdwh8dvgcXtl2GgfKmwA45yi4PSsZP5wxnEdDiILc4Yom/HVbET48UAFH11/oyWnReGzmCNw4OoGnAVPAYBgJQBabHRv2ncVftxWhqLYVAGDQqrBwahoevjaDR0KIQsyZulb8/YsirNtTjk6b89TfkQmR+OGM4bh9UjK0nGGZ/BzDSACpa7Fg/d6z+MeXRagyWwAApjAt7p+WjvuvGYrYSL3kColIpurmDry+owRv7jyDZosNAJBsMuCR64bhO9lDMChCJ7lCokvzWRjZvn07fve736GgoACVlZV4//33MX/+/F6fk5eXh8WLF+Pw4cNITU3Fr371KzzwwAN9fs1gCiPNHVYcPNuEA+VNOFDeiP1lTTjb2O7++mCjAY9cl4GFOWmI1HMaGCI6z9xhxVu7SvHql8Woaba4H0+NCcPElGhkpZgwYUg0JqSY+PeD/EJf998ev1tbW1uRlZWFhx56CN/97ncvu3xxcTFuu+02PProo1i9ejU2b96MRx55BElJSZgzZ46nLx9QOqx2HK4w40B5Iw6WN2F/eSOKaltxqfg3enAUHpqegduzk3mVXSK6JKNBi0dnDMcD1wzF+r1n8cZXJThe1Yyy+naU1bfjowOVAABFAYbHR2JiiglZKdGYmGLCmCQjJ1Qjv3VFwzSKolz2yMgvfvELfPTRRzh06JD7sYULF6KxsREbN27s0+v405ERi82O+tZO1LV0oqGt0/3/+tZO1LV2or7V4v5/aV0bbI6Lf7xDosMwMcXk/iQzPsUEo0ErYW2IKNA1tVtx6Kzzw86BsiYcPNv9aKuLRqUgLSYcMRE6921QhA4x4d3vx3b9G6FT8+wdumI+OzLiqZ07d2L27NndHpszZw6eeuqpHp9jsVhgsZw/BGk2m31S23sF5Th4tgkdVnvXzQGLzflvR9e/FqsdFpsDHVY72q12tHXaPXqNuEgdJnZ9MslKcR4+jWMPCBF5iSlMi+kj4rpNfljTbMHBs85hYOewcCNqWzpRVNvqbo6/HI1KgU6jct7UKmjVKuhd97sec/1fq1ZBrShQqZwfUlWKArUCqBSl677z/66vXxhxXHnH9ej5+xcu45tQxKzV3X9MTsH4ISYpr+3zMHLu3DkkJiZ2eywxMRFmsxnt7e0ICwu76DkrVqzA8uXLfV0a8k7U4N/7Kzx+nkaluD9BuD5RxF7wqSImQo+YCB3SYsORbDLw0wURDaj4KD1uGJ2IG0Y7//YKIVDR1IGy+jY0dB25bWjtRH3X0d361q4jvS3Or1lsDtgcArZOzz+AUeDKThsUvGGkP5YsWYLFixe775vNZqSmpnr9dW4em4j0mHAYtCoYtGroNSrotWoYtGoYNM7HXI87/69CdJgOxjANAwYRBQxFUTAkOgxDoi/+8HcpbZ02NLVbYbUJdNqdR4c7u25Wu/OxTpvD/bjVLuAQAkII2B0CDoGu+85/z98XcF2cWMA5hO1qFHAPaHc9ILrf7RMBvz851K+NTIiU9to+DyODBw9GVVVVt8eqqqpgNBoveVQEAPR6PfR63w9lzMtKxrwsn78MEVFACddpEK7zy8+qFKR8PmPOtGnTsHnz5m6Pbdq0CdOmTfP1SxMREVEA8DiMtLS0oLCwEIWFhQCcp+4WFhaitLQUgHOIZdGiRe7lH330URQVFeHnP/85jh07hj//+c9Yt24dnn76ae+sAREREQU0j8PInj17kJ2djezsbADA4sWLkZ2djaVLlwIAKisr3cEEADIyMvDRRx9h06ZNyMrKwosvvoh//OMfQT/HCBEREfUNp4MnIiIin+jr/ptXWSIiIiKpGEaIiIhIKoYRIiIikophhIiIiKRiGCEiIiKpGEaIiIhIKoYRIiIikophhIiIiKRiGCEiIiKpAuKyjK5JYs1ms+RKiIiIqK9c++3LTfYeEGGkubkZAJCamiq5EiIiIvJUc3MzTCZTj18PiGvTOBwOVFRUICoqCoqieO37ms1mpKamoqysLGiveRPs68j1C3zBvo7Bvn5A8K8j16//hBBobm5GcnIyVKqeO0MC4siISqVCSkqKz76/0WgMyjfYhYJ9Hbl+gS/Y1zHY1w8I/nXk+vVPb0dEXNjASkRERFIxjBAREZFUIR1G9Ho9li1bBr1eL7sUnwn2deT6Bb5gX8dgXz8g+NeR6+d7AdHASkRERMErpI+MEBERkXwMI0RERCQVwwgRERFJxTBCREREUgVdGFm5ciWGDh0Kg8GA3Nxc5Ofn97r8O++8g9GjR8NgMGDChAn4+OOPu31dCIGlS5ciKSkJYWFhmD17Nk6ePOnLVeiVJ+v397//Hddddx0GDRqEQYMGYfbs2Rct/8ADD0BRlG63W265xder0StP1nHVqlUX1W8wGLotE8jbcObMmRetn6IouO2229zL+NM23L59O+bNm4fk5GQoioINGzZc9jl5eXmYPHky9Ho9RowYgVWrVl20jKe/177i6fqtX78eN910E+Lj42E0GjFt2jR8+umn3Zb59a9/fdH2Gz16tA/XoneermNeXt4l36Pnzp3rtlygbsNL/X4pioJx48a5l/GnbbhixQpMnToVUVFRSEhIwPz583H8+PHLPk/2vjCowsjbb7+NxYsXY9myZdi7dy+ysrIwZ84cVFdXX3L5r776CnfddRcefvhh7Nu3D/Pnz8f8+fNx6NAh9zK//e1v8cc//hF/+ctfsGvXLkRERGDOnDno6OgYqNVy83T98vLycNddd2Hr1q3YuXMnUlNTcfPNN+Ps2bPdlrvllltQWVnpvq1Zs2YgVueSPF1HwDlr4IX1nzlzptvXA3kbrl+/vtu6HTp0CGq1GnfccUe35fxlG7a2tiIrKwsrV67s0/LFxcW47bbbMGvWLBQWFuKpp57CI4880m2H3Z/3hK94un7bt2/HTTfdhI8//hgFBQWYNWsW5s2bh3379nVbbty4cd2235dffumL8vvE03V0OX78eLd1SEhIcH8tkLfhH/7wh27rVVZWhpiYmIt+B/1lG27btg2PP/44vv76a2zatAlWqxU333wzWltbe3yOX+wLRRDJyckRjz/+uPu+3W4XycnJYsWKFZdc/s477xS33XZbt8dyc3PFD3/4QyGEEA6HQwwePFj87ne/c3+9sbFR6PV6sWbNGh+sQe88Xb9vstlsIioqSrzxxhvux+6//35x++23e7vUfvN0HV9//XVhMpl6/H7Btg1ffvllERUVJVpaWtyP+ds2dAEg3n///V6X+fnPfy7GjRvX7bEFCxaIOXPmuO9f6c/MV/qyfpcyduxYsXz5cvf9ZcuWiaysLO8V5kV9WcetW7cKAKKhoaHHZYJpG77//vtCURRRUlLifsyft2F1dbUAILZt29bjMv6wLwyaIyOdnZ0oKCjA7Nmz3Y+pVCrMnj0bO3fuvORzdu7c2W15AJgzZ457+eLiYpw7d67bMiaTCbm5uT1+T1/pz/p9U1tbG6xWK2JiYro9npeXh4SEBGRmZuKxxx5DXV2dV2vvq/6uY0tLC9LT05Gamorbb78dhw8fdn8t2Lbhq6++ioULFyIiIqLb4/6yDT11ud9Bb/zM/InD4UBzc/NFv4MnT55EcnIyhg0bhnvuuQelpaWSKuy/SZMmISkpCTfddBN27NjhfjzYtuGrr76K2bNnIz09vdvj/roNm5qaAOCi99yF/GFfGDRhpLa2Fna7HYmJid0eT0xMvGjs0uXcuXO9Lu/615Pv6Sv9Wb9v+sUvfoHk5ORub6hbbrkF//znP7F582Y8//zz2LZtG+bOnQu73e7V+vuiP+uYmZmJ1157DR988AHefPNNOBwOXHPNNSgvLwcQXNswPz8fhw4dwiOPPNLtcX/ahp7q6XfQbDajvb3dK+97f/LCCy+gpaUFd955p/ux3NxcrFq1Chs3bsQrr7yC4uJiXHfddWhubpZYad8lJSXhL3/5C9577z289957SE1NxcyZM7F3714A3vnb5S8qKirwySefXPQ76K/b0OFw4KmnnsL06dMxfvz4Hpfzh31hQFy1l67cc889h7Vr1yIvL69bg+fChQvd/58wYQImTpyI4cOHIy8vDzfeeKOMUj0ybdo0TJs2zX3/mmuuwZgxY/DXv/4Vzz77rMTKvO/VV1/FhAkTkJOT0+3xQN+GoeKtt97C8uXL8cEHH3Trp5g7d677/xMnTkRubi7S09Oxbt06PPzwwzJK9UhmZiYyMzPd96+55hqcPn0aL7/8Mv7v//5PYmXe98YbbyA6Ohrz58/v9ri/bsPHH38chw4dktqD1FdBc2QkLi4OarUaVVVV3R6vqqrC4MGDL/mcwYMH97q8619Pvqev9Gf9XF544QU899xz+OyzzzBx4sRelx02bBji4uJw6tSpK67ZU1eyji5arRbZ2dnu+oNlG7a2tmLt2rV9+sMmcxt6qqffQaPRiLCwMK+8J/zB2rVr8cgjj2DdunUXHQ7/pujoaIwaNSogtl9PcnJy3PUHyzYUQuC1117DfffdB51O1+uy/rANn3jiCXz44YfYunUrUlJSel3WH/aFQRNGdDodrrrqKmzevNn9mMPhwObNm7t9cr7QtGnTui0PAJs2bXIvn5GRgcGDB3dbxmw2Y9euXT1+T1/pz/oBzg7oZ599Fhs3bsSUKVMu+zrl5eWoq6tDUlKSV+r2RH/X8UJ2ux0HDx501x8M2xBwnnZnsVhw7733XvZ1ZG5DT13ud9Ab7wnZ1qxZgwcffBBr1qzpdkp2T1paWnD69OmA2H49KSwsdNcfDNsQcJ6lcurUqT59IJC5DYUQeOKJJ/D+++9jy5YtyMjIuOxz/GJf6JU2WD+xdu1aodfrxapVq8SRI0fED37wAxEdHS3OnTsnhBDivvvuE88884x7+R07dgiNRiNeeOEFcfToUbFs2TKh1WrFwYMH3cs899xzIjo6WnzwwQfiwIED4vbbbxcZGRmivb3d79fvueeeEzqdTrz77ruisrLSfWtubhZCCNHc3Cx++tOfip07d4ri4mLx+eefi8mTJ4uRI0eKjo6OAV+//qzj8uXLxaeffipOnz4tCgoKxMKFC4XBYBCHDx92LxPI29Dl2muvFQsWLLjocX/bhs3NzWLfvn1i3759AoB46aWXxL59+8SZM2eEEEI888wz4r777nMvX1RUJMLDw8XPfvYzcfToUbFy5UqhVqvFxo0b3ctc7mfmz+u3evVqodFoxMqVK7v9DjY2NrqX+clPfiLy8vJEcXGx2LFjh5g9e7aIi4sT1dXVA75+Qni+ji+//LLYsGGDOHnypDh48KB48sknhUqlEp9//rl7mUDehi733nuvyM3NveT39Kdt+NhjjwmTySTy8vK6vefa2trcy/jjvjCowogQQvzv//6vSEtLEzqdTuTk5Iivv/7a/bUZM2aI+++/v9vy69atE6NGjRI6nU6MGzdOfPTRR92+7nA4xH/+53+KxMREodfrxY033iiOHz8+EKtySZ6sX3p6ugBw0W3ZsmVCCCHa2trEzTffLOLj44VWqxXp6eni+9//vpQ/EBfyZB2feuop97KJiYni1ltvFXv37u32/QJ5GwohxLFjxwQA8dlnn130vfxtG7pO8/zmzbVO999/v5gxY8ZFz5k0aZLQ6XRi2LBh4vXXX7/o+/b2MxtInq7fjBkzel1eCOepzElJSUKn04khQ4aIBQsWiFOnTg3sil3A03V8/vnnxfDhw4XBYBAxMTFi5syZYsuWLRd930DdhkI4T2MNCwsTf/vb3y75Pf1pG15q3QB0+73yx32h0lU8ERERkRRB0zNCREREgYlhhIiIiKRiGCEiIiKpGEaIiIhIKoYRIiIikophhIiIiKRiGCEiIiKpGEaIiIhIKoYRIiIikophhIiIiKRiGCEiIiKpGEaIiIhIqv8PdbniMYKQWwoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy                 #loading our favorite library\n",
    "from matplotlib import pyplot    #and the useful plotting library\n",
    "%matplotlib inline\n",
    "\n",
    "nx = 41\n",
    "dx = 2 / (nx - 1)\n",
    "nt = 20    #the number of timesteps we want to calculate\n",
    "nu = 0.3   #the value of viscosity\n",
    "sigma = .2 #sigma is a parameter, we'll learn more about it later\n",
    "dt = sigma * dx**2 / nu #dt is defined using sigma ... more later!\n",
    "\n",
    "\n",
    "u = numpy.ones(nx)      #a numpy array with nx elements all equal to 1.\n",
    "u[int(.5 / dx):int(1 / dx + 1)] = 2  #setting u = 2 between 0.5 and 1 as per our I.C.s\n",
    "\n",
    "un = numpy.ones(nx) #our placeholder array, un, to advance the solution in time\n",
    "\n",
    "for n in range(nt):  #iterate through time\n",
    "    un = u.copy() ##copy the existing values of u into un\n",
    "    for i in range(1, nx - 1):\n",
    "        u[i] = un[i] + nu * dt / dx**2 * (un[i+1] - 2 * un[i] + un[i-1])\n",
    "        \n",
    "pyplot.plot(numpy.linspace(0, 2, nx), u);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Learn More"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For a careful walk-through of the discretization of the diffusion equation with finite differences (and all steps from 1 to 4), watch **Video Lesson 4** by Prof. Barba on YouTube."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/jpeg": "/9j/4AAQSkZJRgABAQAAAQABAAD/2wCEABALDA4MChAODQ4SERATGCgaGBYWGDEjJR0oOjM9PDkz\nODdASFxOQERXRTc4UG1RV19iZ2hnPk1xeXBkeFxlZ2MBERISGBUYLxoaL2NCOEJjY2NjY2NjY2Nj\nY2NjY2NjY2NjY2NjY2NjY2NjY2NjY2NjY2NjY2NjY2NjY2NjY2NjY//AABEIAWgB4AMBIgACEQED\nEQH/xAAaAAEAAwEBAQAAAAAAAAAAAAAAAQIDBAUG/8QAOBAAAgIBAwEFBgUDBAIDAAAAAAECAxEE\nEiExBRNBUXEUIjJhcoEjNJGhsRXB0UJS8PEz4SRik//EABYBAQEBAAAAAAAAAAAAAAAAAAABAv/E\nABcRAQEBAQAAAAAAAAAAAAAAAAABESH/2gAMAwEAAhEDEQA/APvwAAAK7knhtfqBYEEgAQM84AkE\nbl5r9SQAIbS8SHKKaTay+i8wLArGSlna08cPBIEggkAAAAAAAAAAABBJD4WQBJww16d1UZ7UrFJ8\nPmLXOP5OyM4SxtknlZWH1QFgAAIJMtRb3NTnjPK8fmBoDzb+0YR77daod20k0/jyuPTxOzT25pqV\nlkZWSgm2vHjloDcERalFSTynymiQAAAAEAAc1/4mrpqeduHN4fljH8l9PqIajvNieITcG34tDBuQ\nMpdWY13Z1Ftcmko7WufB/wDQG4II3LLWVleAFgZV31WqTrnGSi8PD6FXqqtrkpZSSlleOegGwOO2\n5WQpvjujstSlGXDWeOf1NoylTXZPUTWNzaflHwCNpPCb6/IrVZG2uM4PMZLKEZwnFSjJNNZTMKPw\ntVbSvhaVkflnOV+q/cDqAAUAAAAAAAAAAFLJONcpJZaTaR4uvTrhTOF7bltdiTWGm+rzzjw4PdOO\nfZtE5qTgnj4cpNx9PIC+leO9rzmNcsR+Sx0KPtLTJ4bn/wDnI6Kq41QUI9F4vq/my4GNGqq1Emq3\nLK65i0YXVQXa2nsS9+UJpvPodpSVMZXQted0E0vv/wBAcOl0dF2ihOccTWWpp4cXl8m1Wra7Khqr\nYvd3ak0vFk+wR7pVO23u087U0s/J/I3sqhbTKqS9ySw0vIDzpvVe1J6rupR7ixxjDP8A9eorlfLt\nKKr7pQWmTimnlHWtFF2b52WTfdutZfRP+/BpDTwharFncoKH2QHB2RddHTaSq7Y3Otzbj8sf5Lan\nX2wpUoKEPxZV7pptLDwunmdK0VcI1KuUouqLjF/Jk1aRU1qFdtiSbbbeW8vLyBjPUqqbtnGM5qjc\n5QfXnoiPadTp961Srk9kpxdecceDybLQ0qvu+XHu+75fgS9HCWe8nOeYOvLfRf5Ax0+p1Er61cq+\n7ui5Q2ZyunX9TuMvZ4KVUlnNScY+nH+DUCQAAAAAAADn19Vl+itqqltlOOM/ydBncrHBqqSjLwbW\nQOGWjlfb3+xVOupwqi/Bvxf6L9y2h01lN6bjtrjTGuKznGM5Jzq9232urOcf+Px8uprTDVKxOy+u\ncF1ShhgdQOTXQumodzGcms522bTPRV6iN2ba7Ixx1lapfsB3HJrKrJ3aeyKc41SbcF4vGE/sU1E7\nv6io1Qm9tPH+1tvx9MfuXlXfVp6oQ1CTisSlOO5yYE6XSKt3W2xi7Lp75LHC4wl+iMtbprZ3KdCW\ne7lFfKT6MtBayeduqqeOv4fT9yl8tXRFOzVV89Eqst/uB16WvutNXXjGyKil5JGxjpnKVEZTsjY5\nLKlFYTRsAAAAAAcl25a6CjjdKqai30zlHn6eWo0movhulYlJxjHHDk1uT/do9W+l2OEoS2zg8p4z\nx4omU4wU5SWFFZb8zUrNnXkNXOmNSslZbNb08NKU/wCyWMnaoR9srusqw7IKPK+GS/7/AGN9LK+U\nXK9RW73o48M+BuS0kDzY6WcoapZxfZPLcs4254Xpj+TuV9Tlt3JNycUn4tdS1ltdbSnOMc9MvqRb\n15Sq1N2pdNk3VCcGpKKwmuOhvDSXSjzsi47IpPo9rOydW7UV25+GLXrnH+DGWpjVqLY2ywlt2rxe\nfI1us5jP2Kz2a2Dscp2Pf5JS+XyLSoulVbUpf6ouDlzxw8G89TTW8Tsinz1flyytmsorvjTKWJSW\nenCXzZNq8c1fZrV0LJ2fBBJJeecv/B0xjJ66c8Yiq1HPm8tmEe0q+6tusi41QntUlzuecfyWhrbF\nqY031Rqc1mOJ5F0464yUs4aeOCxzaL4LF5Wz/lnSRoAAAAAAAAAAAAAAAAAAAAACCQBBIAAAAAAA\nAAAAAAABBE5bYSljOFnBYgDx27tRXoVFpTvsVs4xXwR6v98L7nR2NGfcWW2XOx3Tc+VjC6L9kjuh\nVCttwilnyJjGMViKSXyAsAAIOTW2S73T0QS/Fnlya+FR5/XodhSdcbFicU18wPHeonT7fqXaoxbf\ndvHM9sef34+x1XQ1k9Nit1W2t87vd2rHKTR3OuDSThFpdOC0YqKwlhAUoh3dFcGoxcYpYj0XHgaE\nEgAAAAAEHna+UsaxR5xSuG/U9EpKmubm5RT3rbL5osSs9NLUyb7+uuEfDbLJuEklhEkpJjyp0XPU\nNquWKZuyL/3ZeePtk7dRRK9QcbO7a5+FP+TcDTBJqKTeWvHzOW3Qwnf38X+LuTy+emVj9zqAXHBf\n2XDUWynZY8N52pefDX3KWaTOqVSllSpcZSfXblYR6RjKlvV13KWFGMoteecY/gu1nIotDQpZUOM5\n2593PoaWaauxPMcSbT3Lrx0NQNXI8zR6mdfa2p0VsJNNKyFmPdfmvXoemVdcHNTcVuXRliKkAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAgkgAASBAAAkAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAgkAQCQBBIAAAAAAAAAAFLZbK5S8lktF5i\nn5gSAQBIIJAAgkACABIIAEggkAAAAAAAAAAAAAAAAAAAAAAEEkACSCQAAAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAy1H5ez6WXh8EfQpqPy9n0svD4I+gFjl1NlkbIxqznDfCydRXbHduxzjGQMPbIJ\nRzF5fh4on2qDxhPnD/5+hdaepPOznzyV9kow13aw/m/LAGdesXdxdkXuf7v/AIy8tXCGdyaxw388\nZLvTUv8A0IPT1POYJ54Ayjq13soSi854X6f5Fts67nl4io5isfEzWOnqg04wSa+ZMqa5yzKOX6gY\nvWRTw4yclhPHnx/kvK1y07nDcnz0WXwyZ6aubzjHKbw+uC3dQ2KO33V0WQJqblXGUsZay8FyEklh\ncJACQAAAAAAAAAAAAAAAAAAAIAkAAAAAAAAAAQSAAAAAAAAAAAAAAAAAAAAAAAAAAAAGWo/L2fSy\n8Pgj6FNR+Xs+ll4fBH0AsQSeZ2jotRqNVXdTLb3cVj32udyfT0yB6QPL7M0Wr0+tvt1E1Kuz4I72\n9nPTn9TW/sx3doQ1S1VkFHH4a6PH3A7bLYVR3WTjCPm3gd7W7O73x34ztzzg4+1qNRqdL3WmhVJu\nS3d55eOPmZ1aG2Gvjc1BQS3ZT97O3G308QPSBWqU5VxdkVGb6pPoeQuy9Q/aYyn7t8Zp/iSfLeY4\n8gPZJOTQ02U6SNdsVGS44m5fuzPs/s16K2yb1Nl2/wAJeAHY7a1aq3OO9rKjnkV2QsjurkpLplPJ\n52p7Pst7VhqYQhHENrs3Pc1h8Y+5t2XpbdLXNWqCcmsKHThYz9wO4AgCQQSAAAAAAAAAAAAAAQCQ\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABlqPy9n0svD4F6FNR+Xs+ll4fAvQCx5fatus\njOMNPuUXteYxy5PdyvlweoeZ2r2nLRThCuEW3tcpSkkknLH3A20k7vaL4XSlLE245jhKPhyTbZq1\nrIxrrzTxmXBloO1FrdTfQq1F09XuTzz4Fr9fbVro6eOllOLx+InwsgT2m7lVB0uaW579nXo8fvg0\nnK72FxWfae5zlf7sf5Kdp6v2PTKashCTkox3LO5+ReOonPU11RSwob7X5Z6Jfv8AoBxQs1EaapOV\n+O/4TXLhjnP3PUjYpTnFZzDGTy7O151UwtdKkpuxvEktsYvGfma6fteF+ruo7vmnOWpJ558F4gaa\nyesjqao0JOqfEnj4MPl/oYy1GpjptQ47nNW+63D/AEm9utUL9LHcoRuk47ZrEuhrpLnfCzckpQsl\nB48cMBTO6eihOcNtrhlxfmU0M9VOM/aq9jT93pydYA86p6r+qT3Ofdc9fhxhYx885L6qyz2jRyql\nZsc3vSXGMPr98HcAIJAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\ny1H5ez6WXh8EfQpqPy9n0svD4I+gFjOyiq7He1Qnjpuing0Mrb6qMd7ZGG54WX1YEwoqrea6oQfT\nMYpFyN8efeXHXklNNZTygKWVV24VtcJpcrdFPBMa4QbcYpN9cEwnGcd0GmvNDdHONyz5ZAzWloim\nu6g023ys9epaNFUJboVQjLnlRSfJM7IQcVKSTk8L5ss5JdWlgCk6a7JRlOuEpR+FuOWvQmqqFUds\nFhZb+7J3x/3L9SQJAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAAABlqPy9n0svD4F6FNT+Xs+ll4fCvQCxxazQvU6rTXfh4pbeJxy+fI7SAPIt7HsnXdHvofieO18\n+9nL55fgenTUqaI1xUYqKxiKwjyr5ayen7QqjXqlKUvw5JLLXC4/c9HR1uvQ11p2ZUetnxfcCvZ2\nms0lDqsnCS3OS2Rwll5HsFXtntXPeHP2VC7T6ex2wt+JJRk8tvCTfpkjubv6730Xd3bhiSl8C46r\n5gdWsotvVfdWRg4TUnujnOCut0b1VFtWYLe4vLXk0+Rr46lxg9O8xUvxILhyXyZlo4dox1c5amyD\n07ztiuq8gOaXYlruomr4KNM21HZw05Zw/kl0PZPM7Rlr1rtOtKp9xx3uIp8Z8CmuWq7nX93G12Sc\nVUoeWF0++QPWByVSsnr9+2cYOlNxl4PL/ctro6l1RelksqXvR8ZR8Un4MDpB52nr7R9t32Tj7K+k\nH8SLa9WrWaKddds1Gx73Dok01z98Ad4AAkAAAAAAAEAkAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAGWp/LWfSy8fgXoU1P5az6WXj8C9ALEEkAc+tvemqjYkmt8Yy+Sbxkwp107oaOW\nFFX7nL7J8FF2zprYan2Zu2enTc49MpdcP7Guo1+kr1NFF2e8kt8Pdylw+c+mQModoWPRK3MJTV6r\nkl0acscfbk0v1llV9ilFRqhFvfhvHBbQ67R9qVuene9Vy8Y4w/M7AOXs6+eo0yss8W8cc48MnUCQ\nIBIAgEgCASAIJAAAAAAAAAAAAAAAAAAAAAAQBIAAAAAAAAAAAAAAAAAAAAAAAAAAAEASCCQAAAy1\nP5az6WXh8C9Cmp/LWfSy8fhXoBYgkAebPsfs+KtTr2K94niTWec4J9g0He6ex+/JLZXum3nr/wCz\nq1mneopUYy2yjJSi8Z5Tyc9eidMdIlYm6MpuSxuTQDQ6fQ6GrOl2whZLbnPDfTB1u6uMnFzjuSy1\nnk4v6fZ7N3Hexa75WN7cY97c8fc2loYvVPUKct7XR9AOiqyF1cbK5KUJcprxLnPoqHpdNGlz37ej\nxg6AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAgASAAIBJAEgAAAAAAAAAAAAAAAAAAAAAAAy1\nP5az6WXj8K9Cmp/LW/Sy8PhXoBYgk83tDX36XWU1V0boTx776Zz0/TkDzoaHtl23q6+cq7LE04WJ\nYju5x5cGktH2nO/SSsjKaqlFv8ZJLGctrxeDb+papblKuCTk8Sw8RSk1l/bk7ey7bbuz6rL3mx53\nPGPFgdgAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAIBIAAAAAAAAAAAAAAAAAAEASCCQAAA\nAAAAAMtV+Wt+ll4fAvQz1X5W36WaQ+BegFiGk+qJObVa2jSSqjdJqVstsElltgdGF5DhcIxes0yU\n831rZ8XvfD6msJwsgp1yUovo0+GBYHM9dplSrpWxVbnsUn55wbucIw3uSUcZz4AWBSu2u2O6ucZr\nzTyZ1amu2y2Ecp1PEtyx4AbgxlqaIzjCV0FKfMU319BRqadQpdzZGe1tPHg08AbAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAQSAIJAAAAAAAAAAy1X5a36WXh8C9Cmq/LW/Sy8\nPgXoBY5tRpXdfTarXDum3jannJ0kAeZZ2RXKFkFc0pLC91Pas5a+fPmehVDu6owznC64weVLsnUS\nm4ysh3Tkv9csuO/OP04PR01M69HCm2SlJRw2mBzx7NfskqZaiUn3neRm4r3XnPQ6p0Rt0zotbnGU\ndsn0yZaHR+yKa3uW71MnpLv6grlOPdqbk1l55jjGPVIDfR6KjRVuvTxcYt55eSK9LKGovtd8pK3H\nuuKxE6SQOG/s6N9kLJWPdGMY/CvBp/2NdLpPZp2ONjcJyclDC91t5fJ0gAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADHVflbfoZpD4F6Geq/K2/SzSHwr0AsQ\nSQB4ne6yOo1tsbLntmoxi6vdjHzS8TSnVa+zWaaueYRlDNn4XV5458OD0dVqFpoRk1lSmo9cYy8G\nNeu732Vwikr8t5fRJAR2tdbTpVKne5b48QjubWef2N7pWy0kp6ZJ2uOYKXHJzvXyWlVzhHKuVcop\n543bTrus7qmdmM7VnAHP2c9a6Ze3xhGzdxt6YOfS2WS7Q1b/APkRiliKsWYt+a/wdeh1XtdLs2bM\nPGMnNZ2tGG78LMo543Lwko/3Apq9XrK76e5TlXsjKa7ptyy+fTgdn6ntC93q6Ci1/wCLdDapLL5f\nz+RvDXuyFDUMSstdUk30azn+Cs9fKGl1FzjFuiza4p5yuP7MB3s12rVCSt96p7sJ7E/Dn9TvIk9s\nG8N4XReJy9n6162qc3p7ads3HFixnDwB2EHEte5a2zSqn360225cYwsP7/2M32pinSzdWe+Scmnx\nHov5YHog856/UV2ahX6XZXCcYVT3fHl4+xXR9p26iVKlTCKsU23v6bZ7ePMD0wedd2jKtTkoKSc3\nCvLwntWW2/1X2KaTtmN99FEq9s7ao2ZzwsrOAPUJPJ1Has6dPK2NaluUp1pvHux8fub6ftOF+us0\nqhiUFlyzw+nT9QO8EZXHK5AEgAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAQBIIJAAAAAAAAAx1\nX5a36WaQ+Behnq/ytv0M0h8C9ALAEAROEZxcZxUk/BlXRU9i2L8N5jjjD/4y+V5lZThDG6SW54WX\n1YHDXqdDqNJPUVxTqpsefdx7yOy26FVErbXthFZl8kVv00LtNOjG2M14eBedSsolVZzGcdssfNAc\nk+1dDVjfbtzJxXuvqi992mqndvqTdcFOb259P4MV2Lp1CUXbqJblJNynzykn/CN7dFG2Vv4k4q2C\nhLa8PjxyBErdI46dT2R7xqdafHPXJSVum7i2dlO2FdiU014prkvV2dRXXRFqVncfBKby/wBSn9Lr\nVV9attavlunvln9AOx2RUtu5bv8Abnk4v6pRXRO26M6lC3u2ms8v09TssphZltYk1jcuqMNFoKtF\np+5hKycM5/EluArZ2hoouzfNJxT3Zj4Lr/JEdRo7NPVbXGM4Oe2GI45yVu7I0999lsrLouxNOMZ4\nXKSeF9kavQQUMQnPPeq1OTzyAlq6LKdS3FyhQ3GxNeSyyZWaOudVb7uMmnKEcGf9MgnqnG65vU8T\nUpZS9F4cG70tMra7JQTnWsRfkBxT7Q7Lr08YylHuuMLY315X9zoc9LvrrVcX3kHJYj/px1/cwh2J\npoWb1bflPKTnwuGsenLN1oIRdGycl3VbqTzy48ePnwApt0mp0lVm2CqfEFNeXH9g7qYWXruMOmKk\n2kuU/L9CdL2fTpaO6TnbHc5LvXuabKz7Pi777o22qd0NjTlmKXyQHVBwnCEo4ccZiUpthbKzasSh\nLbJfP/ovXBV1xhH4YpJFNPQqXY85lZNyb/58gNgAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAY6v8rb9DNIfAvQz1f5S36GaQ+BegFiCTye1P6p7dp3osdwubPnz0f2ApDsi6Dkpahz\ni5qWXnlKW7n+C39Ks30ylbGXdyjLLT4xnhfqYY7WjGxN3zVjfPu5gtzxt+2Op6PZMbodm0x1MZRt\nSe5TeX1fUDtIJAEAkAQCQBAJAEAkAQCQBAJAEAkAAAAAAAAAAAAAAAgkgASQSABAAkEEgAQAJAAA\ngkgASAAAAAAAACAJAAAAAAABjqvytv0s0j8C9CmpTemtSWW4svH4V6AWAAAAAAAAAAEAACQAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAQCQBBIAEAkAAQAJAAAAAACAJBBIEAkAAAAAAAAAQCQAAAAAAAA\nAAAEAkAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAQCQAAAAAACCSABJBIAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAACCSABJBIAAAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAIJIAZWcZ56kOSXVpFLKd7b3bcxcSi0qWfeys55XTnIG+SHOMcZf\nUw9lb62yfGEadzxFbvdTzjzA0M/aK+/7nd7/AKft6mmDlnoYS7QhrN0lKCxtT4b6Zf2YHUAAJAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAQAAJIJAAgk\nAAAB/9k=\n",
      "text/html": [
       "\n",
       "        <iframe\n",
       "            width=\"400\"\n",
       "            height=\"300\"\n",
       "            src=\"https://www.youtube.com/embed/y2WaK7_iMRI\"\n",
       "            frameborder=\"0\"\n",
       "            allowfullscreen\n",
       "        ></iframe>\n",
       "        "
      ],
      "text/plain": [
       "<IPython.lib.display.YouTubeVideo at 0x7f2744542b38>"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from IPython.display import YouTubeVideo\n",
    "YouTubeVideo('y2WaK7_iMRI')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<link href='http://fonts.googleapis.com/css?family=Fenix' rel='stylesheet' type='text/css'>\n",
       "<link href='http://fonts.googleapis.com/css?family=Alegreya+Sans:100,300,400,500,700,800,900,100italic,300italic,400italic,500italic,700italic,800italic,900italic' rel='stylesheet' type='text/css'>\n",
       "<link href='http://fonts.googleapis.com/css?family=Source+Code+Pro:300,400' rel='stylesheet' type='text/css'>\n",
       "<style>\n",
       "    @font-face {\n",
       "        font-family: \"Computer Modern\";\n",
       "        src: url('http://mirrors.ctan.org/fonts/cm-unicode/fonts/otf/cmunss.otf');\n",
       "    }\n",
       "    div.cell{\n",
       "        width:800px;\n",
       "        margin-left:16% !important;\n",
       "        margin-right:auto;\n",
       "    }\n",
       "    h1 {\n",
       "        font-family: 'Alegreya Sans', sans-serif;\n",
       "    }\n",
       "    h2 {\n",
       "        font-family: 'Fenix', serif;\n",
       "    }\n",
       "    h3{\n",
       "\t\tfont-family: 'Fenix', serif;\n",
       "        margin-top:12px;\n",
       "        margin-bottom: 3px;\n",
       "       }\n",
       "\th4{\n",
       "\t\tfont-family: 'Fenix', serif;\n",
       "       }\n",
       "    h5 {\n",
       "        font-family: 'Alegreya Sans', sans-serif;\n",
       "    }\t   \n",
       "    div.text_cell_render{\n",
       "        font-family: 'Alegreya Sans',Computer Modern, \"Helvetica Neue\", Arial, Helvetica, Geneva, sans-serif;\n",
       "        line-height: 135%;\n",
       "        font-size: 120%;\n",
       "        width:600px;\n",
       "        margin-left:auto;\n",
       "        margin-right:auto;\n",
       "    }\n",
       "    .CodeMirror{\n",
       "            font-family: \"Source Code Pro\";\n",
       "\t\t\tfont-size: 90%;\n",
       "    }\n",
       "/*    .prompt{\n",
       "        display: None;\n",
       "    }*/\n",
       "    .text_cell_render h1 {\n",
       "        font-weight: 200;\n",
       "        font-size: 50pt;\n",
       "\t\tline-height: 100%;\n",
       "        color:#CD2305;\n",
       "        margin-bottom: 0.5em;\n",
       "        margin-top: 0.5em;\n",
       "        display: block;\n",
       "    }\t\n",
       "    .text_cell_render h5 {\n",
       "        font-weight: 300;\n",
       "        font-size: 16pt;\n",
       "        color: #CD2305;\n",
       "        font-style: italic;\n",
       "        margin-bottom: .5em;\n",
       "        margin-top: 0.5em;\n",
       "        display: block;\n",
       "    }\n",
       "    \n",
       "    .warning{\n",
       "        color: rgb( 240, 20, 20 )\n",
       "        }  \n",
       "</style>\n",
       "<script>\n",
       "    MathJax.Hub.Config({\n",
       "                        TeX: {\n",
       "                           extensions: [\"AMSmath.js\"]\n",
       "                           },\n",
       "                tex2jax: {\n",
       "                    inlineMath: [ ['$','$'], [\"\\\\(\",\"\\\\)\"] ],\n",
       "                    displayMath: [ ['$$','$$'], [\"\\\\[\",\"\\\\]\"] ]\n",
       "                },\n",
       "                displayAlign: 'center', // Change this to 'center' to center equations.\n",
       "                \"HTML-CSS\": {\n",
       "                    styles: {'.MathJax_Display': {\"margin\": 4}}\n",
       "                }\n",
       "        });\n",
       "</script>\n"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from IPython.core.display import HTML\n",
    "def css_styling():\n",
    "    styles = open(\"../styles/custom.css\", \"r\").read()\n",
    "    return HTML(styles)\n",
    "css_styling()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "> (The cell above executes the style for this notebook.)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
